function [maxn, poly, midres, trares, simpres]=question4(fid)
% QUESTION4 - fid 为输出文件的文件句柄。
%   

    n = 16;
    N = [];
    i =1;
    resp =1;
    while 1
        N =[N,n];
        mid = multimid(@fun,0,1,n);
        res =abs(mid-pi);
        midres(i) =res;
        if res > resp
            maxn(1) = n;
            break;
        end
        resp = midres(i);
        i= i+1;
        n = 2*n;
    end
    l = length(N);
    N = 1./N;
    x = log(N(1:l-2));
    y = log(midres(1:l-2));
    poly(1,:)=polyfit(x,y,1);
    
    n = 16;
    N = [];
    i =1;
    resp =1;
    while 1
        N =[N,n];
        mid = multitra(@fun,0,1,n);
        res =abs(mid-pi);
        trares(i) =res;
        if res > resp
            maxn(2) = n;
            break;
        end
        resp = trares(i);
        i= i+1;
        n = 2*n;
    end
    l = length(N);
    N = 1./N;
    x = log(N(1:l-2));
    y = log(trares(1:l-2));
    poly(2,:)=polyfit(x,y,1);
    
    n = 16;
    N = [];
    i =1;
    resp =1;
    while 1
        N =[N,n];
        mid = multisimp(@fun,0,1,n);
        res =abs(mid-pi);
        simpres(i) =res;
        if res > resp
            maxn(3) = n;
            break;
        end
        resp = simpres(i);
        i= i+1;
        n = 2*n;
    end
    l = length(N);
    N = 1./N;
    x = log(N(1:l-2));
    y = log(simpres(1:l-2));
    poly(3,:)=polyfit(x,y,1);
    

    n = 16;
    N = [];
    i =1;
    resp =1;
    while 1
        N =[N,n];
        mid = Romberg(@fun,0,1,n,2);
        res =abs(mid-pi);
        rom2res(i) =res;
        if res > resp
            maxn(4) = n;
            break;
        end
        resp = rom2res(i);
        i= i+1;
        n = 2*n;
    end
    l = length(N);
    N = 1./N;
    x = log(N(1:l-2));
    y = log(rom2res(1:l-2));
    poly(4,:)=polyfit(x,y,1);
    
    n = 4;
    N = [];
    i =1;
    resp =1;
    while 1
        N =[N,n];
        mid = Romberg(@fun,0,1,n,3);
        res =abs(mid-pi);
        rom3res(i) =res;
        if res > resp
            maxn(5) = n;
            break;
        end
        resp = rom3res(i);
        i= i+1;
        n = 2*n;
    end
    l = length(N);
    N = 1./N;
    x = log(N(1:l-2));
    y = log(rom3res(1:l-2));
    poly(5,:)=polyfit(x,y,1);
    
    n = 4;
    N = [];
    i =1;
    resp =1;
    while 1
        N =[N,n];
        mid = Romberg(@fun,0,1,n,4);
        res =abs(mid-pi);
        rom4res(i) =res;
        if res > resp
            maxn(6) = n;
            break;
        end
        resp = rom4res(i);
        i= i+1;
        n = 2*n;
    end
    l = length(N);
    N = 1./N;
    x = log(N(1:l-2));
    y = log(rom4res(1:l-2));
    poly(6,:)=polyfit(x,y,1);
    
    adpt = adp(@fun,0,1);
    adpres = pi-adpt;
  
fprintf(fid,'P107第四题:\n');
fprintf(fid,'复合中点积分公式的误差(n=16,32,...)\n');
fprintf(fid,'%12.6e &',midres);
fprintf(fid,'\n\n');
fprintf(fid,'复合梯形积分公式的误差(n=16,32,...)\n');
fprintf(fid,'%12.6e &',trares);
fprintf(fid,'\n\n');
fprintf(fid,'复合simpson积分公式的误差(n=16,32,...)\n');
fprintf(fid,'%12.6e &',simpres);
fprintf(fid,'\n\n');

fprintf(fid,'Romberg积分公式的误差(k=2,n=16,32,...)\n');
fprintf(fid,'%12.6e &',rom2res);
fprintf(fid,'\n\n');
fprintf(fid,'Romberg积分公式的误差(k=3,n=4,8,...)\n');
fprintf(fid,'%12.6e &',rom3res);
fprintf(fid,'\n\n');
fprintf(fid,'Romberg积分公式的误差(k=4,n=4,8,...)\n');
fprintf(fid,'%12.6e &',rom4res);
fprintf(fid,'\n\n');
fprintf(fid, ...
 ...
'当n达到下面数时因舍入误差，精度不再增加。（依次是中点，梯形，simpson,Romberg 2,3,4）\n');
fprintf(fid,' %d &',maxn);
fprintf(fid,'\n\n');

fprintf(fid,'各个求积公式的误差对h的双对数线性拟合结果:\n');
fprintf(fid,'%12.8e &%12.8e \n',poly');
fprintf(fid,'\n\n');
fprintf(fid,'%12.8e  ',exp(poly(:,2)));
    
fprintf(fid, '\n\n自适应求积法的结果:');
fprintf(fid,'%4.20f ',adpt);
fprintf(fid,'\n\n Matlab 给的近似值%4.20f',pi);

    
end

function f = fun(x)
% FUN - 
%   
    sx = size(x);
    f = ones(sx);
    
    for i = 1:sx(1)
        for j = 1:sx(2)
            f(i,j)=4/(1+x(i,j)^2);
        end
    end
    
end



